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pj Abstract 

^ Almost all sources of high energy particles and photons are associated with jet phenomena. 

I Prominent sources of such highly relativistic outflows are pulsar winds, Active Galactic Nuclei 

H and Gamma-Ray Bursts. The current understanding of these jets assumes diluted plasmas which 

^ are best described as kinetic phenomena. In this kinetic description particle acceleration to 

I ^1 ultra-relativistic speeds can occur in completely unmagnetized and neutral plasmas through 

insetting eff"ects of instabilities. Even though the morphology and nature of particle spectra are 

^ understood to a certain extent, the composition of the jets is not known yet. While Poynting- 

1^ flux dominated jets (e.g. occuring in pulsar winds) are certainly composed of electron-positron 

QQ plasmas, the understanding of the governing physics in AGN jets is mostly unclear. 

^ In this article we investigate how the constituting elements of an electron-positron-proton 

^ plasmas behave differently under the variation of the fundamental mass-ratio mp/nte. 

(^ We studied initially unmagnetized counterstreaming plasmas using fully relativistic three- 

^ dimensional particle-in-cell simulations to investigate the influence of the mass-ratio on particle 

►>■ acceleration and magnetic field generation in electron-positron-proton plasmas. We covered a 

S^ range of mass-ratios mp/nie between 1 and 100 with a particle number composition of np+/ne+ 

of 1 in one stream, therfore called the pair-proton stream. Protons are injected in the other one, 

therfore from now on called proton stream, whereas electrons are present in both to guarantee 

charge neutrality in the simulation box. We find that with increasing proton mass the instability 

takes longer to develop and for mass-ratios > 20 the particles seem to be accelerated in two 

phases which can be accounted to the individual instabilities of the different species. This means 

that for high mass ratios the coupling between electrons/positrons and the heavier protons, 

which occurs in low mass-ratios, disappears. 

Subject headings: methods: numerical; acceleration of particles; galaxies: jets; plasmas; insta- 
bilities 



1. Introduction 

Radiation observed from astrophysical systems like GRBs or AGN usually possesses a nonthermal 
emission spectrum. This is believed to arise from particle acceleration in the vicinity of relativistic shocks 
or within the counterstreaming plasma itself. 

In recent PiC simulations it has been shown that most particle acceleration occurs within the jet 
dNishikawa er¥l[2003| [2005] [2006] [2008| [20091 [Chang et al.|[2008| [Spitkovsky[[2008al [Dieckmann eTaL 
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Ramirez -Ruiz et al. 2007 Silva et al. 2003 1 and is mostly caused by plasma instabilities like the Weibel 
(Weibel 1959 1 or twostream instability ( [Buneman 1958| ). Both instabilities create current filaments with 
surrounding magnetic fields and are therefore a plausible source for particle acceleration and the generation 
of observed long-lasting magnetic fields. Particle acceleration can also occur along with shocks where first 
order Fermi acceleration ( [Fermi[|1949[ ) is assumed to be the relevant process, which was shown by kinetic 
simulations only recently ( Spitkovsky|2008b I. 



In the present work we focus on the main properties of the plasma instabilities and describe the in- 
fluence of the fundamental mass-ratio trip/ nig in mixed electron-positron-proton compositions by means of 
relativistic three-dimensional simulations of counterstreaming plasmas. 

The paper is organized as follows: In chapter [2] the underlying code is described briefly, in section [3] 
we illustrate the setup of the performed simulations. In chapter |4] we present the results of our simulations 
which we discuss and draw some conclusion in section [5] 



2. Description of the code 

Particle-in-Cell simulations are an essential tool in understanding relativistic coUisionless plasma physics. 
Therefore we developed a three-dimensional fully relativistic MPI-parallelised PiC code called ACRONYM 
(Another Code for moving Relativistic Objects, Now with Yee lattice and Macroparticles). Maxwells equa- 
tions are evolved in time by employing a second-order leapfrog scheme (see e.g. Taflove and Hagness|2005| l. 
The particles affect the electromagnetic fields through charge currents which are deposited on the grid by us- 
ing a second-order Triangular Shaped Cloud (TSC) scheme (see e.g. [Hockney and Eastwood) 1988"! ) adopted 
from [Esirkepov ( [2001 1. The particles are moved via a second-order force interpolation within the Boris 
push ( Boris|1970 ). In order to guarantee that the divergence of the magnetic field remains close to zero, the 
electric and magnetic fields are stored in the form of a staggered grid, the so-called Yee-lattice ( ,Yee,1966^ . 
With this setup the code is second-order both in space and time. 

Extensive tests of the code have been successfully completed from which we conclude that the total 
relative error in energy conservation is less that 3 x 10"^ and the divergence of the magnetic field stays below 
a value VB/B < 10 /Ad in the simulated space for all times. 



3. Simulation setup 

In the simulations presented here we use two counterstreaming plasma populations, one representing 
the background medium consisting of 6e^ and 6p''' per cell (proton stream) and the other incorporating the 
jet containing 4e~, le'^ and 2/?^ per cell (pair-proton stream), from which we find the background density 
ratio Ujet/rihg = 2/3 and the ratio np+/ne+ = 1 in the pair-proton stream. In the lab frame (the rest frame of 
the simulation box) the two streams are counterstreaming along the z-direction with a Lorentz factor y = 10 
(yS = v/c = 0.995) each, the electron distribution has a thermal velocity of Vth,e - 0.1c in every direction in 
the restframe of the moving medium, the thermal velocity of the protons is Vth,p = 0.1c • {mp/me)~^. This 
setup resembles situations as they are believed to exist in jets running into the interstellar or intergalactic 
medium. 

Three-dimensional simulations with five different compositions of counterstreaming plasmas using 
128 X 128 X 512 cells with a total of 167 million particles (20 particles per cell) and mass-ratios nip/me 
between 1 and 100 have been performed. In addition to that another simulation with a mass-ratio of 100 
with twice the number of cells in each of the perpendicular directions (and therefore four times more par- 
ticles) has been performed in order to show the influence of the periodic boundary conditions. As pointed 
out by Fonseca et al. ( |2008 l), 20 particles per unit cell (on average) in combination with quadratic parti- 



cle interpolation are sufficient to eliminate most of the numerical noise. Nevertheless, we have performed 
simulations up to twice the numbers of particles per cell and no significant changes were observed. 

Periodic boundary conditions have been applied in all three dimensions. Due to the quick development 
of the selfconsistent electromagnetic fields it is redundant to solve Poissons equation at the initial time and 
the fields can be initialized with zero without loss of generality. 

The cell size is set to be equal to the Debye-length of the plasma. Ax = Ajj - (X-qT I Ajirii.e^)^!'^ , and 
the timestep is restricted by the CFL-criterion, c • Af < Ax/ V3. This results in a At between 0.035w~^ 
and 0.050(jjp^ (normalized to the plasma frequency Up = {Ane^nlnie)'^^^) and the cellsize Ax ranges from 
2x10^ cm to 3 X 10"^ cm depending on the different mass-ratios employed. The simulations were evolved 
for 2500 to 4500 timesteps to roughly 120 to 220 ojZ} (the exact numbers for each simulation can be found 
in Table [T]l. The characteristic scales of interest in a counterstreaming electron-positron-proton plasma are 
of the order of several proton skin depths, cjcopi - {ympC^/Ane^ny^^, for a plasma with the density n and 
the average proton energy ynipC^. 



4. Simulation results 

The results of our simulations can be divided into two main findings: (1) new insights into the evolu- 
tionary behavior of twostream instabilities in multi-component plasma and (2) the change in the distribution 
function of the particles, in particular the acceleration caused by the instability. 

Due to the huge amount of data, simulation results have been written every tenth step for the fields 
(electric and magnetic fields, currents) and every hundredth step for particle data. 



4.1. Analysis of the fields and currents 

In this section we analyze the evolutionary behaviour of the electric and magnetic fields and the currents 
in the simulations conducted. From previous simulations of filamentation instabilities in pair plasmas (see 



e.g. Silva et al. 2003 1 it is well known how magnetic and electric fields evolve. We compare the behaviour of 



plasmas with different mass-ratios. The most significant quantity in this context is the transverse magnetic 
field energy averaged over the entire computational domain B\ = {B^ + B^), since strong magnetic fields are 
essential to create and maintain the flux tubes observed in kinetic instabilities, furthermore the point in time 
the instability peak occurs and also the existence of a second peak, respectively. 

In Fig. [T] we therefore compare the time evolution of the transverse magnetic field energy B^ computed 
in the lab frame as function of different mass-ratios nip I nig . 

It is evident that the maximum value of the transverse magnetic field energy reached in the different 
simulations are comparable, even though it has to be noted that for the non-pair plasma the maximum energy 
can only be found in the second peak. The development of the second peak shows a nicely observable 
dependence on the fundamental mass ratio: For the lower mass-ratios no two peak structure can be seen, 
while with increasing mass-ratio a clear distinction can be made. 

When looking at the time until the instability fully develops, one can see that for higher mass-ratios 
mp/nif, it takes longer to reach the peak value for the magnetic field. If one compares simulations with mass- 
ratios of 1 and 100 one can explain what is happening in the counterstreaming plasma: First an electron- 
positron instability develops almost simultaneously for both mass-ratios (single peak for mass-ratio 1 and 
first peak for mass-ratio 100). If a third and heavier species exists another peak will be apparent at later 
times (which can be seen in Fig. [T]l. This behavior is not observable for medium mass-ratio simulations 
since both peaks overlap and can not be distinguished anymore. 

The existence of two instabilities in the plasma has important impact on the amplitude and duration 
of the instability: Clearly the instability lasts longer for the high mass-ratio simulations, since in this case 
the heavy protons are accelerated slower compared to the lighter electrons/positrons but are able to stabilize 
the flux tubes for a longer period of time. Another result to note is that the maximum amplitude decreases 
with increasing mass-ratio. This effect can be attributed to the lower number of particles constituting each 
instability. 



mp/nie 


At[cOp'] 


Ax [cm] 


number of cells 


number of ion skin depths 


1.0 


0.0498 


23225 


128x128x512 


120 X 120 X 480 


5.0 


0.0385 


23107 


128x128x512 


57 X 57 X 228 


20.0 


0.0360 


21614 


128x128x512 


28x28x112 


42.8 


0.0355 


21338 


128x128x512 


20 X 20 X 80 


100.0 


0.0408 


21198 


128x128x512 


13x13x52 


100.0 


0.0574 


29830 


256x256x512 


36x36x72 



Table 1: Setup of the simulations performed 
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Fig. 1 . — Comparison of the time evolution of the transverse magnetic field energy for different mass-ratios. 
With increasing mass-ratio the instability takes longer to develop. For the two highest mass-ratios one 
notices that the instability develops in two phases. 



To inspect the effect of employed mass-ratios on the temporal evolution of instabilities, we looked at 
the nature of the flux tubes more closely. The flux tubes in the simulations with mass-ratios of mp/nie = 5.0 
and mp/rrie = 100.0 are illustrated in Fig.s|2]and|3j respectively, which show the particle number density 
in particles per cell at three different locations perpendicular to the direction of streaming. The pictures in 
Fig. |2] are chosen such that the uppermost row shows the onset of the instability and the lowest pictures 
are roughly taken at the time the maximum of the instability occurs (cf. Fig. [T]). In Fig. |3]the upper row 
of slices is taken at the moment the electron/positron instability peaks, the second set of pictures show 
the time between the two instabilities (compare with the black curve in Fig. [llat 80oj"') and the last set 
shows the point in time when the proton instability reaches its peak. Both simulations show the archetypical 
behavior of filamentation instabilities: Flux tubes develop, which in turn merge until only two flux tubes 
survive. But for the high mass-ratio simulations this whole process happens twice. In an early stage (which 
resembles the first peak of the instability in Fig. [TJ flux tubes arise. In a later stage (second peak) flux tubes 
of different strengths exist, one of them is more pronounced. The explanation is that during the first stage of 
the instability the flux tubes are carrying more of the lighter particles. The second stage is then associated 
with a flux tube of heavier particles which takes longer to develop, but is also able to exist for a much longer 



timespan. 

The combination of Fig. s[T] [2] and [3] suggests that the instability is evolving in two phases. In the first 
phase light particles are accelerated and in the second phase the heavier particles are also involved in the 
instability. For mass-ratios of nip/nie < 20.0 the instabilities and therefore particle acceleration of the two 
species cannot be separated. There are two possible reasons for this: Either the coupling of the two species 
is still strong enough to co-accelerate the heavier particles or the time scales of the two instabilities are still 
matching rather well. A strong argument for the second option is the fact that the instability time scale 
for the heavy species increases less than hnear with the mass-ratio for a certain size of the computational 
domain. 

Considering only mass-ratios mp/nie > 20.0 the instabilities are clearly separated which is also apparent 
in the double-hump structure in the energy-diagrams of Fig. [T] and in Fig. [3] illustrating the development of 
the flux tubes in two phases. 

In the larger simulation (twice the size in the perpendicular directions) it takes even longer for the 
proton instability to develop because the flux tubes have more space to develop and therefore it takes longer 
for them to merge until only two are left. This is also the reason for the slightly different slopes in the two 
simulations with mass-ratio 100. The maxima and minima of the magnetic energy occur around the point in 
time, when the current density in the direction of streaming averaged over the whole computational domain 
changes its sign. When (in the larger simulation) the proton instability kicks in (at around 65a)~^) the flux 
tubes are not yet fully merged down to two and therefore the proton instability does not grow with the same 
rate as in the smaller simulation. The resulting flux tubes around the maximum of the second instability still 
resemble the two flux tube regime as seen in the smaller simulations. 



4.2. Particle distribution 

As described in chapter |3] all simulations were initialized with thermal particle distributions (width of 
the thermal distribution 0.1c and 0.1c • (nip/nig)"^ for electrons and protons, respectively) which are then 
boosted with a Lorentz factor of y = 10 in either direction. While some particles gain a lot of energy during 
the simulation in total, the shape of the particle distribution is also changing. To analyze and quantify the 
change of the particle energy distributions we utilize two distinct types of graphs: (1) a two-dimensional 
plot in the lab frame relating the absolute value of the momentum parallel (v||/c • y with vy = \v^\) and 
perpendicular (v^/c • y with Vx = (v^ -i- Vy)^) to the initial streaming direction, respectively and (2) an 
one-dimensional plot of the distribution of the particles speed in the lab frame. 

In Fig.|4]we show the time evolution of the electron and positron distribution (all electrons and positrons 
in both streams are plotted in a 2D histogram for the mass-ratios 5 (upper panels) and 100 (lower panels). 
The same conjuncture is illustrated in Fig.|5]for the protons, but the axis are different. 

As expected, the particle distribution in the early stage of the simulation (before the onset of the first 
instability hump) is centered at the initial Lorentz boost y = 10 with a thermal width of 0.1c for the elec- 



trons/positrons and 0.1c • (mp/me)~^ for the protons. Several electrons have already been accelerated and are 
streaming towards higher transverse velocities. At later times in the low mass-ratio simulation one can see 
that both electrons and protons are getting accelerated simultaneously, as illustrated in the next two images 
in the upper rows of Figs. |4] and [5] As emphasized before there is no clear distinction between a low and 
high mass instability in this case, thus the simultaneous acceleration is clearly in agreement with the energy 
evolution (see Fig.[T]l. 

When going to higher mass-ratios and early times (lower leftmost image of both figures) only the 
electrons are visible, as the protons still remain at their initial momentum. When the initial phase of the 
instability is over {t k 10co~^) the electrons are no longer significantly accelerated. The protons are still 
gaining more energy until about t = 1 10a>~^ which is roughly at the maximum of the proton instability (see 
Fig.[T]for comparison). This supports the idea of two almost separated instabilities. 

In both cases the particle distributions show a diffusion-like behaviour: The initial distribution in par- 
allel direction is stretched from y w 5 to y w 15. The particle energy is converted into perpendicular field 
energy which in turn accelerates the particles. We want to stress an interesting feature of the different par- 
ticle distributions in the diverse simulations: Both simulations have a proton distribution which is strongly 
elongated in the perpendicular direction and an electron distribution whose center is below the initial y = 10 
and extends more or less equal in the parallel and perpendicular direction. In the high mass case it is obvious 
that the electrons are partially decelerated to lower energies. 

In Fig. [6] we show a one-dimensional plot of the temporal evolution of the total momentum distribution 
of all the particles (from both proton and pair-proton stream) in the lab frame. Particles moving in negative 
(positive) z-direction are plotted with a negative (positive) total momentum and the narrow peak for each 
time shows the protons, the broad ones represent the electrons/positrons. One can see that most of the 
acceleration of the electrons and positrons is happening between lOco'^ and 70a>~^ the protons still gain 
more energy until about 1 10a;„^ . This corresponds with the behaviour of the particles seen in Fig. Hand in 
Fig-lU 



5. Discussion and summary 

We have conducted several simulations of counterstreaming plasmas with different mass-ratios 1 .0 < 
nip/nie < 100.0 in order to investigate the influence of the mass-ratio on the development of twostream 
instabilities. We can draw two major conclusions from our simulations: (1) the physics of acceleration in 
mixed counterstreaming plasmas can be understood by using particle-in-cell simulations and (2) we are able 
to show that there is indeed a strong implication of the mass-ratio on the results but in order to extrapolate 
the behaviour to the physical mass-ratio one still needs to perform some simulations with higher mass-ratio 
of 200 or 500. 

We are able to demonstrate that the mass-ratio has a qualitative and not only quantitative effect on 
the simulated physics. For very low mass-ratios only one instability develops which changes for mass-ratios 
around 50. Taking these findings as a starting point we conclude that further simulations of counterstreaming 



plasma have to be conducted with mass-ratios of 100 and above to find the correct physical behaviour. A 
quantitative change has also been found: For mass-ratios > 20 the instability is developing more slowly but 
also lasts longer (cf. Fig. [TJ because the flux tubes can be sustained for a greater timespan. This can be 
explained by a two stage instability, that can be seen especially in the highest mass-ratios. Here in a initial 
phase a pair instability develops as in pair-only plasmas resulting in several flux tubes. In some of the flux 
tubes the (lighter) electrons and positrons are streaming and these ones are developing much stronger. After 
reaching its maximum the light-particle instability is decreasing and some other flux tubes (carrying mostly 
the heavier protons) are growing stronger (see Fig.|3]). 

Besides the investigation of the nature of the instability itself the acceleration of particles was a subject 
of research. From two-dimensional plots of particle momentum parallel versus perpendicular to the original 
streaming direction in the lab rest frame we concluded that most of the particle acceleration happens in 
the transverse direction. Furthermore, the two stage process can also be observed here: While light and 
heavy particles are accelerated almost simultaneously for the low mass-ratio simulations, the electrons are 
accelerated stronger and earlier compared to the protons in the high mass-ratio simulations. Additionally, 
at a later stage the electron distribution stays almost the same while the much heavier protons still gain 
energy (see Fig.s|4]and|5]l. This gives important insight into the fundamental acceleration mechanism: Only 
particles involved in the instability may be accelerated. 

Since the instability itself and the restrictions on numerical parameters (i.e. the mass-ratio) are now 
better understood, it is necessary to conceive the implications on the underlying physics of astrophysical jet 
phenomena. From this study it is yet not possible to impose strong limits on the possible composition of jets, 
since the major difference between pair and mixed plasma seems to be the timescale of how the instability 
develops. Obviously an extensive parameter study is necessary to cover the full range of possible physics. 
Therefore in a next stage the influence of the composition of the background and jet plasma on the particle 
acceleration shall be investigated. 

TB thanks the Deutsches Zentrum fiir Luft- und Raumfahrt and the LISA-Germany collaboration for 
funding. OF and UG are grateful for funding from the Elite Network of Bavaria and FS would like to thank 
the Deutsche Forschungsgemeinschaft (DFG SP 1124/1). Simulations have been conducted at Hochleis- 
tungsrechenzentrum Stuttgart under grant "iotmrofi". 
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Fig. 2. — Development and merging of the flux tubes for the simulation with mplrtie = 5.0 and a resolution 
of 128 X 128 X 512 cells. The colors pertain to the particle density (in particles per cell) shown at three 
different locations along the the direction of streaming (at 100, 300 and 500 Ajc). One can see the fluxtubes 
developing and merging over time. Both fluxtubes have about the same strength, which shows that the 
instabilities of the electrons/positrons and the protons are still srongly coupled for this mass-ratio. 



12 



t = 40(o'' 



density [particles/cell] 
t = 40 (o'^ 




100 




100 



10 



t = 40 to" 




100 



I 



10 



t = 80(o'' 



t = 80 ffl' 



t = 80(o'' 




100 



10 




100 



10 




100 



10 



t = 200(o"' 



t = 200 (o'^ 




100 



10 




100 



i 



10 



t = 200(o'' 




100 



10 



Fig. 3. — Development and merging of the flux tubes for the simulation with Hip /me = 100.0 and a resolution 
of 256 X 256 x 512 cells. The colors pertain to the particle density (in particles per cell) shown at three 
diff'erent locations along the the direction of streaming (at 100, 300 and 500 Ax). In the upper three pictures 
it is to be seen that first the flux tubes develop, later they nearly vanish (which corresponds to the dip in 
the magnetic field energy in Fig. Ml in the black curve at SOcoZ,^) and then the flux tube in the lower right 
corner grows stronger (lower set of pictures). This can be attributed to the independent instabilities at high 
mass-ratios. 
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Fig. 4. — Development of the absolute value of the momentum distributions (v^/c • y over vy/c • y) of all the 
electrons and positrons parallel and perpendicular to the direction of the original flow in the lab rest frame. 
First one can state that most of the particle acceleration is happening in the transverse direction. Secondly, 
for a mass-ratio of m^/m^ = 5.0 electrons and positrons are accelerated until the simulation stops, while for 
the higher mass-ratio {mp/nie = 100.0) the electrons stop gaining energy around about t = lOco'^. 
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Fig. 5. — Development of the absolute value of the momentum distributions (v^/c • y over vy/c • y) of all the 
protons parallel and perpendicular to the direction of the original flow in the lab rest frame. First one can 
state that most of the particle acceleration is happening in the transverse direction. Secondly, for a mass-ratio 
of nip /me = 5.0 the protons are getting accelerated over the entire timespan, while for the higher mass-ratio 
{nip/ Me = 100.0) their acceleration starts at a later time during the simulation. 
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Fig. 6. — Time evolution of all the particles (both proton and pair-proton stream particles) momentum 
distribution in the lab frame for a mass-ratio of mp/rrie = 100.0. Particles moving to the left (proton stream) 
are shown with a negative total momentum, particles moving to the right (pair-proton stream) with a positive 
one. The narrow peak for each time shows the protons, the broad ones represents the electrons/positrons. 



